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FOR INVISCID DUCT FLOW PROBLEMS 



A. Verhoff 



McDonnell Aircraft Company 
St. Louis, Missouri 



March 1990 



ABSTRACT 



Highly accurate far field computational boundary conditions for inviscid^ 
two-dimensional isentropic duct flow problems are developed from analytic 
solutions of the linearized, second-order Euler equations. The Euler 
equations are linearized about a constant pressure, rectilinear flow 
condition. The boundary procedure can be used with any numerical Euler 
solution method and allows computational boundaries to be located extremely 
close to the nonlinear region of interest. Numerical results are presented 
which show that the boundary conditions and far field analytic solutions 
provide a smooth transition across a computational boundary to the true far 
field conditions at infinity. The cost of upgrading first-order boundary 
conditions to second-order is slight. 
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INTRODUCTION 



Numerical solution procedures for nonlinear fluid dynamic equations 
usually use one or more artificial computational boundaries located at some 
distance from the primary region of interest in order to limit the physical 
domain to finite size. If the flow crossing such a boundary (either inflow or 
outflow) is subsonic, then some type of computational boundary conditions must 
be imposed which simulate the influence of the true far field conditions at 
infinity. These boundary conditions must be such that waves crossing the 
boundary do not produce erroneous reflections back into the computational 
field to degrade the calculations. It is generally acknowledged that simply 
imposing free stream conditions (or conditions at infinity) at computational 
boundaries is usually inappropriate because of the spurious reflections back 
into the computational domain which result. Standard practice has consisted 
of locating the boundaries quite far from the region of interest in an attempt 
to simplify the boundary condition models and minimize any effects of incon- 
sistent modeling. The net effect is a significant increase in the number of 
grid points required for an accurate flowfield calculation. 

A boundary modeling procedure for two-dimensional internal flows was 
presented in Reference 1 which alleviates the difficulties mentioned above and 
also allows the computational boundaries to be located much closer to the 
nonlinear region of interest. The procedure is limited to steady, inviscid 
flow, although the flow can be rotational. It represents a logical first- 
order extension of the so-called characteristic (or zero-order) boundary 
conditions commonly used with inviscid numerical solution methods. Extension 
to axisymmetric or three-dimensional flows is straightforward. 

The analysis presented here extends the first-order analysis to second 
order for isentropic duct flow and provides a logical extension, in an 
asymptotic sense, of the first-order analysis. It also illustrates a consis- 
tent procedure for coupling linearized analytic solutions with nonlinear 
numerical solutions by means of computational boundary conditions. The 
greater accuracy of the second-order boundary procedure allows the computa- 
tional boundaries to be placed even closer to the nonlinear computational 
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region thereby further reducing the number of grid points needed for the 
numerical solution. Even though the analysis and boundary conditions derived 
from it are limited to isentropic flow, many non-isentropic duct flow problems 
have isentropic flow crossing the inflow boundary at which these second-order 
boundary conditions can be imposed. 

The present analysis is based on the Riemann variable formulation of the 
Euler equations given in Reference 2 simplified by the isentropic assumption. 
This represents a natural starting point because the characteristic (or 
zero-order) boundary conditions mentioned above are expressed in terms of 
Riemann variables. The equations are linearized about a constant pressure, 
rectilinear flow condition, which truly represents conditions at infinity. 
These linearized equations are assumed applicable in the far field region 
beyond a computational boundary. Within the nonlinear computational domain, 
strong entropy-producing (i.e., rotational) effects can exist which create 
variations in density, velocity, etc. in the downstream far field in the 
streamline-normal direction which are not necessarily small perturbations. 

Such variations were modeled in the analysis presented in Reference 1. 

However, because of the two-level perturbation procedure and subsequent 
approximate solution method used in that analysis, inclusion of entropy 
effects in the second-order analysis does not seem justified. As mentioned 
above, however, the second-order isentropic analysis is still applicable to 
non-isentropic problems having isentropic inflow conditions. 

The first-order linearized equations analyzed in Reference 1 are homo- 
geneous and were solved using separation of variables and Fourier analysis. 

The second-order equations treated herein are non-homogeneous, but can also be 
solved by the same techniques. These equations are written in different forms 
for the upstream and downstream far field regions because the downstream flow 
is driven by the distribution of flow angle on the computational boundary 
while the upstream flow is driven by the distribution of Riemann variable 
associated with upstream propagating waves. The assumption is made that the 
far field flow is confined by parallel walls, which is not restrictive for 
duct flow problems. 
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The first- and second-order solutions are coupled to the nonlinear 
numerical solution to provide a smooth transition across the computational 
boundary to the true far field conditions at infinity. The coupling is 
accomplished by the boundary conditions. The first-order boundary conditions 
provide distributions of flow quantities to be imposed along the boundary, not 
constant conditions. The second-order boundary conditions provide high 
accuracy corrections to the first-order distributions. The first- and 
second-order boundary conditions represent logical asymptotic extensions of 
the zero-order (or characteristic) conditions. Furthermore, the boundary 
analysis can be coupled with any inviscid numerical solution method. 



Numerical results are presented for a duct/nozzle configuration. Results 
obtained using the second-ot der boundary conditions are compared with those 
using zero- and first-order conditions. An extremely fine grid was used for 
the numerical solution in order to accurately assess the improvement produced 
by the second-order analysis. It was found that the computational boundaries 
could be located very close when the second-order boundary conditions were 
used with no loss in numerical solution accuracy. The reduced size of the 
computational field further reduced the number of grid points needed for the 
numerical solution. The additional computational effort required to upgrade 
the boundary conditions to second-order is very slight. The benefits of the 
improved boundary conditions are also achieved when a coarser grid is used. 
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H . SECOND-ORDER PERTURBATION EQUATIONS 



The system of two-dimensional, steady, linearized Euler equations which 
describe second-order perturbations from a constant pressure state are derived 
in this section. For reasons outlined in Reference 1, the Riemann variable 
formulation of Reference 2 will again be used because of its close 
relationship with characteristic (or zero-order) boundary conditions commonly 
used in numerical solution of the nonlinear Euler equations. As explained in 
the previous section, the analysis will be limited to isentropic flow in ducts 
having parallel walls. 

The two-dimensional, isentropic form of the Euler equations is 
(Reference 2) 



ag 

3t 



+ (q + a) 



^ = 
as 




( 1 ) 



at 
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(q - a) II = + qa 
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Velocity magnitude and speed of sound are denoted by q and a, respectively, 
and P is the logarithm of pressure. The Riemann variables Q and R are defined 
as 



Q = q 3 

(4) 

R = q - a 
y-1 

The flow angle is 0, time is denoted by t, and local distances along and 
normal to the streamline direction are denoted by s and n, respectively. 
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For steady flow the analysis can be greatly simplified by defining a new 
dependent variable 

T E Q - R (5) 

and replacing equations (1) and (2) by 



(m 2 - 1) |I + 2 q M = 0 (6) 

a2 + ^ q2 = 1 (7) 



The local Mach number is denoted by M. Equation (6) is obtained by 
subtracting equations (1) and (2). Equation (7) is obtained by adding 
equations (1) and (2) and integrating. The constant of integration, which is 
proportional to stagnation temperature, can be set to unity by proper choice 
of non-dimensional izing quantities. The simplified form of the steady Euler 
equations is then 



(«' - » II * ^ Q H = 0 



( 8 ) 



m2 



^ _2_ 191= 

9S Y“1 T 9n 



(9) 



a2 + q2 = 1 



( 10 ) 



In regions of the flowfield where nonlinear effects are weak, the flow 
can be treated as a perturbation to a constant pressure, rectilinear flow. 
Such regions occur near and beyond far field computational boundaries. The 
dependent variables in equations (8) and (9) can then be expanded in 
asymptotic series 



T = T„ + Ti + T2 + ... 

e = e* + 01 + 02 + • • • 



( 11 ) 
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The flow direction at infinity is assumed constant and denoted by 0a,; T„ is 
also constant because the flow is assumed isentropic. The perturbation 
quantities T^- and 0 ^ vanish at infinity. For duct flow between parallel 
walls, the far field flow angle 0a, can be assumed zero. 

Consistent with the expansions (11), spatial derivatives in equations (8) 
and (9) can be approximated by 



where x and y are reference Cartesian coordinates and 0i and 02 are measured 
from the x axis which is aligned with the axis of the duct (see Figure 1). 

If expansions (11) and (12) are introduced into equations (8) and (9), 
the resulting first-order perturbed Euler equations are 




( 12 ) 




(13) 




(14) 



The second-order equations are 




(15) 





( 16 ) 
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Velocity and Mach number at infinity are constant and denoted by q„ and M„, 
respectively. The parameter b is defined by 



6 i n/i-M„2 

Procedures for calculating far field quantities such as q„ and M„ are outlined 
in Reference 1. 

Asymptotic expansions of the Riemann variables Q and R can also be 
defined as 

Q = Q® + Ql + Qz + • • • 

(18) 

R = R„ + Rj + R2 + ... 

Using the definition (5) and the expansion (11) for T, it follows that 

T<o = Qo, ~ Rod 

Ti = Qi - Ri (19) 

Tz = Qz - Rz 



Introducing the expansions (11) and (18) into the algebraic equation (10) 
gives the first- and second-order relationships 

(1+H„) Qi - (1-M.) Ri = 0 (20) 

(1+M„) QZ - R2 + 0 + ^ m!) Rf = 0 (21) 

These relations will be used later in Section IV where the boundary conditions 
are derived. 
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Equations (20) and (21) can also be used to express the perturbation 
equations in different dependent variables. One such form, which will be 
convenient later in Sections III and IV, uses Rj and R 2 in place of Tj and T 2 . 
In terms of these variables, the first-order equations are 



(1-M„) ^ + q„ — = 0 



301 



, V 901 9Rl 

Qco (1+M„) ^ ^ = 0 



The second-order equations are 



( 22 ) 

(23) 



(i-M„) ^ " I ^ 



(1+Mco)qo 



(24) 



Qcn (1+M„) 



m! 2 



- i n r— 0 " 

3x ■ ay 2 3y •■1-M„ 

- (1 + 2M„ + ^ M^) fif] 

(1+M„)2q£ 



(25) 



The second-order equations using either choice of dependent variables consist 
of a coupled pair of non-homogeneous equations whose right hand sides (RHS) 
are functions of the first-order solution variables. 
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III. Sii .UTION OF SECOND-ORDER EQUATIONS 



Solutions of the second-order perturbation equations for duct flows are 
developed in this section. The first-order duct flow analysis of Reference 1 
showed that the behavior of the flow in the downstream far field is governed 
by the distribution of flow angle e on the downstream computational boundary. 
This distribution is obtained from the nonlinear numerical solution. Also, 
the first-order downstream boundary conditions for the computational domain 
were expressed in terms of the perturbation variable Tj. Therefore, the 
second-order equations (15) and (16) are appropriate for analyzing the down- 
stream region. 

In the upstream far field region, the flow behavior is governed by the 
distribution of the Riemann variable R on the upstream computational boundary. 
This distribution is obtained from the nonlinear numerical solution. The 
first-order upstream boundary conditions for the computational domain were 
expressed in terms of the perturbation variables 9i and Q^. Therefore, the 
second-order perturbation equations (21), (24) and (25) are appropriate for 
analyzing the upstream region. 



Downstream Region 

The second-order equations (15) and (16) can be combined into the single 
equation 



9^02 9^02 

9x2 9y2 



92 

9x9y 



[mEg? + (1 

4q£ 



+ ^ m!) T?] 



(26) 



From Reference 1, the solution of the first-order system (13) and (14) for 
duct flows can be expressed as 
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( 27 ) 



Gj = I Ap Sin nny 
1 



Tl 




E Ap COS niry 
1 



where Ap represents the Fourier coefficients of the distribution of flow angle 
0 on the downstream computational boundary located at x=0. The width of the 
duct has been taken as unity. Note that the solution (27) for Tj does not 
include a zero-mode term (i.e., Aq) which is related to the mean value of the 
perturbation T-T„ at the boundary. It was pointed out in Reference 1 that 
such a term is a second-order effect. This will be verified below. Using the 
first-order solution in the RHS, equation (26) becomes 



o 

o 9^05 9^0p a><» v-T? O') 

— o- + — o- = 5 E E An, An [(2 + ^ M„) (m^-n^) sin(m-n)iiy 

3yZ 23-^ 



9x^ 



1 1 



(28) 



+ m£ (m+n)2 sin(m+n)Tiy] 



To satisfy the duct wall boundary conditions of zero flow angle, a 
particular solution of equation (28) must have the form 



02 



26^ 



E E Hmn(x) sin(m±n)iry 
1 1 



(29) 



This results in the two component ordinary differential equations (depending 
on choice of sign) for each value of m and n 



- T’^(m+n)2 H^n = ^ m£ (m+n)2 a„, Ap e ("^+0)^^/B 
b 2 Hmn - ^2(m-n)2 = (2 + ^ M^) (m^-n^) a„, A„ (31) 
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The respective solutions of these equations are 




-(m+n)nx/B 



(32) 



(33) 



Referring to the first-order solution (27), the solution for ej can be 
used to exactly represent the computational boundary distribution of flow 
angle obtained from the nonlinear numerical solution. This requires that 02 
vanish at x=0. This can be enforced by adding a solution of the homogeneous 
portion of equation (31) to the solution (33). The full second-order solution 
for 0 in the downstream region is then 



The solution for T 2 can be obtained from equation (15) using the solution 
(29) for 02 which was determined above. After integrating, the full 
second-order solution for T in the downstream region is 




-(m+n)nx/B 



(34) 



+ k2 E E 
1 1 



m^-n^ 

mn 




where 



■ 8 





(35) 
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T = T„ 



2qa>Mo> ” 

— r — E Ap cos nny e 
1 



+ K] E E An, Ap cos(m+n)ny [1 + ^ (m+n)x] 
1 1 ^ 



( 36 ) 



* >:2 E E ^ ^ C0S(m-n),y [(m-n)2 e-'"*"’’*/® - |m2-n2| 

limn 

+ E E An, An [K 3 cos(m-n)ny + K 4 cos(m+n)iry] 

1 1 



The coefficients are defined by 



Kl 

<2 

K 3 

K 4 



]f±i 




4 g4 



1 

4 34 



(2 + 



2 





(1 + 






yzl 

4 




Poo^oj 



(1+ 




+ 



T±i 

4 




(37) 



The solution (36) contains longitudinally decaying plane waves 
corresponding to m=n. They are represented by 



T 2 = K 3 E e 
1 



-2nTTx/6 



qa.Mg, 




g-2nnx/B 



(38) 



At x=0 the quantity T 2 represents the mean value of the perturbation T-Tco at 
the boundary. 
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Upstream Region 



The second-order equations (24) and (25) can be combined into the single 
equation 



From Reference 1, the solution of the first-order system (22) and (23) can be 
expressed as 



where Cp represents the Fourier coefficients of the distribution of the 
Riemann variable perturbation R-R„ on the upstream computational boundary 
located at x=0. As explained in Reference 1, the mean value of R-R<„ on the 
boundary represented by the zero-mode coefficient Cg is a second-order effect 
and is not included in the first-order solution (40). Using the first-order 
solution in the RHS, equation (39) becomes 




(39) 




(1 + 2M„ + ^ m2) Ri]) 



Rl = E Cp COS niiy e 
1 



OtiX/B 



(40) 





46^ Qod(1+Mcd) 



[(1 + M<„ + m£ + m2 + m2) 



r E C„, Cf, (m+n)2 cos(m-n)ny 




1 1 



+ (y+1) m2 e e Cm Cp (m+n)2 cos(m+n)Tiy] 



1 1 
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Following the solution procedure used for the downstream region, a 
particular solution of equation (41) must have the form 



“ 462g.O*H.) ? 5 



( 42 ) 



This results in the two component ordinary differential equations (depending 
on choice of sign) for each value of m and n 

^mn “ Ti^(m+n)^ = (y+ 1) (m+n)^ Cjn Cp ( 43 ) 

Gmn - n^(m-n)2 G^n = [bi(m+n)2 + b2(m-n)2] Cn, Cp (44) 

The coefficients b^ and b^ are defined by 



bj = 1 + Moo + m£ + m2 + m2 

b2 = - (1 + M„ - 3 m£ + m2 - ^ m2) 

The respective solutions of these equations are 

J 

_ 1+1 [If? / ^ p P (m+n)irx/B 

Gmn = (m+n)'^ + b 2 (m-n)'^] “jjj e 



(45) 



(46) 

(47) 



Referring to the first-order solution (40), with the exception of the 
mean value, the solution for Rj can be used to exactly represent the 
computational boundary distribution of R-R„ obtained from the nonlinear 
numerical solution. This requires that all Fourier modes of R 2 except the 
zero mode vanish at x=0. This can be enforced by adding a solution of the 
homogeneous portion of equation (44). The solution for R 2 in the upstream 
region is then 
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R 2 = <5 X E E (m+n) C^) Cp cos(m+n)Tiy e 
I 1 



(mtn)«/6 ^ ^2n.x/6 

1 



g|m-n|nx/6j 

g|m-n|nx/6j 



where 

<5 

<6 

K? 



(y+1 )tiM„ 

863q„(l+M„) 

Tiet (1 + ^ m2) 

(1 * H. - 3H^ * ^ ^ h!) 



'^6 ” cos(m-n).y - 






E E 
1 1 
m?!n 

CO oo 

E E 
1 1 



(m-n)^ 



mn 



Cm Cn cos(m-n)Tiy [e 



(m+n)nx/B 



( 48 ) 



(49) 



The corresponding solution for Q 2 can be obtained from the second-order 
relation (21) as 



Q2 



1-M„ 

1+M„ 



R2 



M„ 



qo.(l+M„)- 



(1 + ^ m 2) fif 



(50) 



The solution for 
(48) for R 2 which was 
second-order solution 



02 can be obtained from equation (24) using the solution 
determined above. After integrating, the full 
for 0 in the upstream region is 



IF 
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0 = 



Z Cn sin nny + IC 3 E E Cn sin(m+n)7iy 



QojC l+Moj) J j J 



k 4 X E E (m+n) C„, Cp sin(nfi+n)7iy 



1 1 



( 51 ) 



- k 5 E E !!!^ C, Cn sin(.-n).y [(-n.n)l^ - (m-n) 

m/n 



k6 n M c„ c„ sin(m-n).y [(m*n) - |m-n| el"-"!"*'®] 



1 1 



The coefficients are defined by 



k3 E ^ (1 + Mo,) + ^ m2] 

4Bq£(l+Mco)2 



k4 = 



TlMa 



l±i _ 

B^q2(1+Moo)^ 



ks = — ^ — (1 + m£ + ^ m2) 

16Bq£(l+Mo,) 



(52) 



ke = 



16Bq2(l+M„)2 



(1 + M,„ - 3m2 + ^ m2 - ^ m2 ) 



The solution (48) also contains longitudinally decaying plane waves 
corresponding to m=n. They are represented by 



D n V V ^2nTiX/B 

K2 = 4 K6 E Cn e 



(53) 



— ^ (1 + m2 + ^ m 2) E Cn 

4B^qcn ^ 1 
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From equatiofi (50) the corresponding plane waves for Q are 



QZ = 



4q„(l+M„)3 



(1 - - 



Id. 

2 



id 

2 



h1) 



,2 „2nnx/B 



(54) 



At x=0 the quantities Q 2 and R 2 represent the mean values of the perturbations 
Q-Q„ and R-Roo at the boundary, respectively. Note that the mean value of the 

R perturbation is predicted by the analysis and not determined directly from 
the numerical solution as are the other Fourier modes. 
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IV. BOUNDARY CONDITIONS DEVELOPMENT 



Second-order boundary conditions for duct flows are developed in this 
section based on the linearized Euler solutions obtained in the previous 
section. They represent a logical extension of the first-order boundary 
conditions developed in Reference 1. 

Downstream Boundary 

The downstream boundary conditions are derived from the definition (5) 
rewritten as 

R = Q - T (55) 

Using the solution (36) for T and applying the relation (55) at the boundary 
(assumed located at x= 0 ), the distribution of R on the boundary (i.e., the 
boundary conditions) is calculated according to 

a 2 gcoMco 

Rb = Qnum " ~V~ J 

CO CO 

+ (K] + K 4 ) EE An, Ap cos(m+n)ny (56) 

1 1 

, n [K2 feiOi-r l!!!L!i ? l , K3I A„ A„ cos(m-n).y 
1 1 

The boundary distribution of Q obtained from the nonlinear numerical solution 
is denoted by Qpum- zero Fourier mode corresponding to m=n in- this 

expression is proportional to the coefficient K 3 . It represents the mean 
value of the R perturbation at the boundary and is clearly a second-order 
effect. 
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Upstream Boundary 



The upstream boundary conditions are derived by applying the perturbation 
relation (50) and the solutions (48) and (51) at the boundary (assumed located 
at x=0). From these relations the distributions of Q and 9 on the boundary 
(i.e., the boundary conditions) are calculated according to 



Qb 



\ CO CO n 

+ TTm^ ^n cos nny + 4 K5 E Cp] 



— 3 (1 + ^ m£) [E Cn COS ni.y]2 



(57) 



and 



0b 



q ( ^ + 1^3 ^ ^ Cn, Cn sin(m+n)ny 



^ -(m-n)] Cm Cn sin(m-n)i,y 



-f kc E E ^ „ 

^ j j mn m-n 



(58) 



+ ke T T ^ [m+n - |m-n|] C^ Cn sin(m-n)ny 



1 1 



The zero Fourier mode proportional to K 5 in equation (57) is again evident. 
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V. APPLICATIONS 



Numerical solutions of the Euler equations have been calculated for 
two-dimensional, isentropic, steady duct/nozzle flow using the second-order 
boundary condition procedures developed in the previous section. Second-order 
boundary condition results are compared with those produced using the conven- 
tional zero-order characteristic boundary conditions and also the first-order 
boundary conditions developed in Reference 1. 

A solution algorithm was used for the nonlinear Euler equations (1) - (3) 
which is based on the method presented in Reference 2. It uses explicit time 
integration to relax to steady state conditions. It should be emphasized that 
the boundary condition analysis is independent of the choice of inviscid, 
nonlinear solution method. 

The duct/nozzle geometry is shown schematically in Figure 1. The flow is 
characterized by Pa,, the downstream pressure at infinity, which produces a 
mass flow per unit area w through the duct. The linearized solutions given by 
equations ( 34 ), ( 36 ), ( 48 ), (50) and (51) are assumed valid in the semi- 
infinite regions I and III and the computational boundary conditions are 
applied at the upstream and downstream boundaries AA and BB of the nonlinear 
computational region II. 

The actual shape of the duct/nozzle and the computational grid are shown 
in Figure 2. The shape is identical to that used for the calculations pre- 
sented in Reference 1 . The nozzle contour is sinusoidal and symmetric about 
the centerline. The computational grid for this portion of the nozzle had 
dimensions 81 x 41 , which is twice the size of that used in Reference 1 . The 
increased grid dimensions allowed second-order effects to be quantified more 
accurately. The area ratio of the nozzle is .75 and the upstream and down- 
stream areas are equal. For these constant area sections of the duct, addi- 
tional rectangular grid cells could be added without altering the basic 81 x 
41 grid. This served to minimize the effect of grid changes on the calcula- 
tions when the computational boundaries were moved in order to assess the 
accuracy of the boundary conditions. 
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Because the configuration is symmetric, calculations were limited to the 
lower half of the duct and a centerline symmetry condition was used. Although 
the configuration used for these calculations is simple, the boundary condi- 
tion analysis of the previous section is general and applicable to unsymmetric 
configurations having unequal upstream and downstream areas. Use of the 
simple configuration is sufficient to demonstrate the validity of the 
analysis. 

The second-order upstream and downstream boundary conditions are given by 
equations (56), (57) and (58). The associated analytic far field solutions 
are given by equations (21), (34) and (36) for the downstream region and by 
equations (48), (50) and (51) for the upsteam region. The zero-order (or 
characteristic) boundary conditions consist of imposing the constant value of 
Q„ and a zero value of 9 along the upstream boundary and the constant value of 
R„ along the downstream boundary. The first-order upstream boundary condi- 
tions consist of imposing additional distributions of Qj and 9i along the 
boundary as determined from equations (20) and (40). The first-order down- 
stream boundary conditions consist of imposing an additional distribution of 
Rl along the boundary as determined from equations (27) and (55). Additional 
details are given in Reference 1. 

Computational results are presented for a single value of Pa,, but with 
the computational boundaries located at two different longitudinal stations. 
One boundary is sufficiently far removed so that all three sets of boundary 
conditions produced essentially the same results within the computational 
domain. The second boundary is extremely close to the nozzle portion of the 
duct so that the relative accuracy of the various boundary conditions can then 
be evaluated. 

Results obtained using the complete grid shown in Figure 2 are presented 
in Figures 3 and 4. This grid has 40 columns of grid cells in both the 
upstream and downstream constant area portions of the duct. As mentioned 
above the results were nearly identical for all three sets of boundary condi- 
tions. Figure 3 shows pressure and Mach number distributions along the 
centerline and lower wall of the duct/nozzle. Pressure, Mach number, and flow 
angle contours are presented in Figure 4. These results serve as a reference 
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for evaluating the accuracy of solutions where the computational boundaries 
are moved closer to the nozzle portion of the duct. 

Results for the shortened computational domain obtained using the zero- 
order boundary conditions are presented in Figures 5 and 6. There was only 
one column of grid cells in both the upstream and downstream constant area 
portions of the duct for this case. Pressure and Mach number distributions 
are shown in Figure 5 and contours are shown in Figure 6. The degradation in 
the solution as a result of the boundary proximity is very evident. 

Application of the first-order boundary conditions to the short computa- 
tional domain gives the results shown in Figures 7 and 8. Eleven Fourier 
modes were needed to accurately describe the boundary distributions of 0 and R 
for this extremely close boundary location. Pressure and Mach number 
distributions are shown in Figure 7 and contours in Figure 8. Linearized 
first-order solution results obtained from equations (20), (27) and (40) have 
been added upstream and downstream of the computational boundaries. 

Significant improvement is provided by the first-order boundary conditions. 
Some mismatch at the boundary is evident, but this can be attributed to 
neglect of second-order interactions. 

Finally, results for the shortened computational domain using the 
second-order boundary conditions are presented in Figures 9 and 10. Eleven 
Fourier modes were also used in obtaining these results. Pressure and Mach 
number distributions are shown in Figure 9 and contours in Figure 10. The 
results within the numerical solution portion of the domain are nearly 
identical to those shown in Figures 3 and 4. Linearized second-order solution 
results obtained from equations (48), (50) and (51) and from equations (21), 
(34) and (36) have been added upstream and downstream of the computational 
boundaries, respectively. It is evident that the linearized far field 
analytic solutions provide for a smooth transition across the computational 
boundary to the true far field conditions at infinity. 

A more quantitative comparison between the three sets of boundary condi- 
tions can be obtained by examining the distribution of flow variables along 
the computational boundaries. Figure 11 shows a comparison of the distribu- 
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tion of pressure and flow < ngle along the upstream and downstream boundaries 
for the short computational domain. Results at the same longitudinal loca- 
tions taken from the reference numerical solution using the long grid (Figure 
3) are also shown. These results exhibit classical asymptotic convergence to 
the reference solution. The larger increment between the zero-order and 
first-order pressure results at the upstream boundary is primarily due to the 
fact that the zero-order boundary conditions impose a zero flow angle along 
the boundary. 

It was pointed out in Section III and also in Reference 1 that the mean 
values of the perturbations at the computational boundaries are neglected 
within the first-order analysis. In Section III they were shown to appear as 
second-order effects described by equations (38) and (53) evaluated at x-0. 

In the downstream region equation (38) indicates that the relationship 



oo 



[E A 
1 



2 

n 




(T-T„) dy = - 



q»M® 

6« 



(1 ♦ ’'f 



(59) 



should hold. Likewise, from equation (53) the corresponding relation for the 
upsteam region is 




(R-Ro.) dy = 



I 

4e2q„ 



(1 + h£ + ^ m£) 



(60) 



These relations are shown in Figure 12 for P® = 0.90. The Fourier 
coefficients and integrals are determined at various longitudinal upstream and 
downstream locations using flowfield data taken from the reference long grid 
solution of Figure 3. Distance measured from the start of the constant area 
portions of the duct is denoted by These results confirm the second-order 
nature of the perturbation mean values. 
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It can be seen from equations (56), (57) and (58) that very little 
additional computational effort is required to upgrade the boundary conditions 
from first-order to second-order. Most of the computational effort at the 
boundaries is spent determining the Fourier coefficients which are needed to 
impose first-order boundary conditions. Therefore, computational boundaries 
having isentropic inflow can be upgraded easily even though a source of 
entropy production may exist within the nonlinear computational domain. 

The results presented in Figures 3 through 12 were obtained using a very 
fine computational grid in order to accurately quantify the differences 
between the boundary condition models. To demonstrate that the improvement 
provided by the second-order procedure can also be achieved with a moderately 
sized grid, the calculations were repeated with the grid used in Reference 1. 
One column of grid cells was used in both the upstream and downstream constant 
area portions of the duct. These results are shown in Figures 13 and 14. 
Pressure and Mach number distributions are shown in Figure 13. Distributions 
of pressure and flow angle along the upstream and downstream computational 
boundaries are compared in Figure 14. Note that these boundaries are located 
slightly upstream and downstream, respectively, from those of the fine grid 
solution. Asymptotic convergence to the reference solution (Reference 1) is 
again evident. Eleven Fourier modes were needed to resolve the boundary 
distributions of 0 and R. 
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VI. 



SUMMARY 



Highly accurate second-order far field boundary conditions for 
two-dimensional, isentropic duct flow have been developed from analytic 
solutions of the second-order linearized Euler equations. These boundary 
conditions represent a logical extension, in an asymptotic sense, of the 
first-order boundary conditions derived in Reference 1 from analytic solutions 
of the first-order linearized Euler equations. The first-order boundary 
conditions of Reference 1 are themselves a logical extension of the zero-order 
(or characteristic) boundary conditions commonly used in numerical solutions 
of the nonlinear fluid dynamic equations. The boundary conditions and 
analytic solutions provide a smooth transition across a computational boundary 
to the true far field conditions at infinity. The boundary procedure is 
general in that it can be used in conjunction with any numerical solution 
method. 

The second-order linearized Euler equations are non-homogeneous and have 
been solved exactly using separation of variables and Fourier analysis. This 
procedure was used in Reference 1 to solve the linearized, homogeneous 
first-order equations. Extension of the second-order analysis to allow 
non-isentropic flow conditions is probably not justified because solution of 
the non-isentropic equations requires an approximate two-level perturbation 
procedure (see Reference 1 for details). However, even if the flow being 
analyzed has an entropy source in the nonlinear region, the second-order 
boundary conditions can still be applied at the inflow computational boundary 
if the incoming flow is isentropic. 

Use of zero-order (or characteristic) boundary conditions requires that 
the computational boundaries be located far from the nonlinear region of the 
flow. Closer placement of the boundaries may result in a significant amount 
of solution degradation. The first-order boundary conditions derived in 
Reference 1 allow the boundaries to be located much closer thereby reducing 
the number of grid points needed for the numerical solution and also the 
number of iterations for solution convergence. This leads to a significant 
reduction in the amount of computational effort required for the nonlinear 
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numerical solution because the additional calculations required for the 
first-order boundary conditions is modest. 

Use of the second-order boundary conditions allows the computational 
boundary to be placed even closer with a further reduction in the number of 
grid points. The amount of additional computational and implementation effort 
is very slight so that isentropic inflow boundary conditions can be upgraded 
with little penalty. 

In Reference 1 it was pointed out that the mean values of the flow 
variable perturbations at the computational boundary were not described within 
the first-order analysis. It was further postulated that these mean values 
represented a second-order effect. The nature of this second-order interac- 
tion was clarified by the present analysis and verified by numerical results. 
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Figure 1. Duct/Nozzle Schematic 




Figure 2. Duct/Nozzle Geometry and Grid 
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Figure 3. Pressure and Mach Number Distributions 
Zero-, First- and Second-Order Boundary Conditions 
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Figure 4. Pressure, Mach Number and Flow Angle Contours 
Zero-, First- and Second-Order Boundary Conditions 
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Figure 6. Pressure, Mach Number and Flow Angie Contours 
ZerO’Order Boundary Conditions 
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Figure 7. Pressure and Mach Number Distributions 
First-Order Boundary Conditions 
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Figure 8. Pressure, Mach Number and Flow Angle Contours 
First-Order Boundary Conditions 
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Figure 9. Pressure and Mach Number Distributions 
Second-Order Boundary Conditions 
p = 0.9 
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Figure 10. Pressure, Mach Number and Flow Angle Contours 
Second-Order Boundary Conditions 
p = 0.9 
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Figure 11. Pressure and Flow Angle Distributions 
on Computational Boundaries 
P„ = 0-9 
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Figure 12. Comparison of Numerical and Theoretical 
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Figure 13. Pressure and Mach Number Distributions 
Second-Order Boundary Conditions 
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Figure 14. Pressure and Flow Angle Distributions 
on Computational Boundaries 
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